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Abstract 

SHARE is a collection of programs designed for the statistical anal- 
q ■ ysis of particle production in relativistic heavy-ion collisions. With the 

. physical input of intensive statistical parameters, it generates the ra- 

tios of particle abundances. The program includes cascade decays of 
all confirmed resonances from the Particle Data Tables. The complete 
treatment of these resonances has been known to be a crucial factor 



behind the success of the statistical approach. An optional feature im- 
plemented is a Breit-Wigner type distribution for strong resonances. 
An interface for fitting the parameters of the model to the experimen- 
tal data is provided. 
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PROGRAM SUMMARY 



Title of the program: SHARE, 



July 2004, version 1.1 



Computer: 

PC, Pentium III, 512MB RAM 



not hardware dependent 



Operating system: 

Linux: RedHat 6.1, 7.2, FEDORA etc. 



not system dependent 



Programming language: F0RTRAN77: g77, f 77 

as well as Mathematica, with temperature-fixed particle yields and excluding 
option of finite width; 

Size of the package: 645 KB directory including example programs (82KB 

compressed distribution archive), without libraries (see 

ht t p : / / wwwasdo c . web . cern . ch / wwwasdo c/minuit / minmain . ht ml 

ht t p : / / wwwasd . web . cern . ch/ wwwasd / cernlib . ht ml 

for details on library requirements ). 

Distribution format: tar gzip file 

Number of lines in distributed program, including test data, etc: 14893 

Keywords: relativistic heavy-ion collisions, particle production, statistical 
models, decays of resonances 

Computer: Any computer with an f77 compiler 

Nature of the physical problem: 

Statistical analysis of particle production in relativistic heavy-ion collisions 
involves the formation and the subsequent decays of a large number of res- 
onances. With the physical input of thermal parameters, such as the tem- 
perature and fugacities, and considering cascading decays, along with weak 
interaction feed-down corrections, the observed hadron abundances are ob- 
tained. SHARE incorporates diverse physical approaches, with a flexibility 
of choice of the details of the statistical hadronization model, including the 
selection of a chemical (non)equilibrium condition. SHARE also offers eval- 
uation of the extensive properties of the source of particles, such as energy, 
entropy, baryon number, strangeness, as well as the determination of the 
best intensive input parameters fitting a set of experimental yields. This al- 
lows exploration of a proposed physical hypothesis about hadron production 
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mechanisms and the determination of the properties of their source. 



Method of solving the problem: 

Distributions at freeze-out of both the stable particles and the hadronic res- 
onances are set according to a statistical prescription, technically calculated 
via a series of Bessel functions, using CERN library programs. We also have 
the option of including finite particle widths of the resonances. While this is 
computationally expensive, it is necessary to fully implement the essence of 
the strong interaction dynamics within the statistical hadronization picture. 
In fact, including finite width has a considerable effect when modeling directly 
detectable short-lived resonances (A(1520), K*, etc.), and is noticeable in fits 
to experimentally measured yields of stable particles. After production, all 
hadronic resonances decay. Resonance decays are accomplished by addition 
of the parent abundances to the daughter, normalized by the branching ra- 
tio. Weak interaction decays receive a special treatment, where we introduce 
daughter particle acceptance factors for both strongly interacting decay prod- 
ucts. An interface for fitting to experimental particle ratios of the statistical 
model parameters with help of MINUIT [1] is provided. 

The x 2 function is defined in the standard way. For an investigated 
quantity / and experimental error Af, 

^2 (/experiment /'theory) 

(A /statistical + A/ S yst cmatic) 2 
NboF A^} a t a points A^f rcc parameters • (2) 

(note that systematic and statistical errors are independent, since the sys- 
tematic error is not a random variable). 

Aside of x 2 the program also calculates the statistical significance [2] , de- 
fined as the probability that, given a "true" theory and a statistical (Gaus- 
sian) experimental error, the fitted x 2 arises at or above the considered value. 
In the case that the best fit has statistical significance significantly below 
unity, the model under consideration is very likely inappropriate. In the 
limit of many degrees of freedom (Ndof), the statistical significance function 
depends only on x 2 /N DoF , with 90% statistical significance at xV^doF ~ 1) 
and falling steeply at x 2 /N Do f > 1- However, the degrees of freedom in fits 
involving ratios are generally not sufficient to reach the asymptotic limit. 
Hence, statistical significance depends strongly on x 2 an d Nd q f separately. 
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In particular, if N DoF < 20, often for a fit to have an acceptable statistical 
significance, a xV^doF significantly less than 1 is required. 

The fit routine does not always find the true lowest y 2 minimum. Specifi- 
cally, multi-parameter fits with too few degrees of freedom generally exhibit a 
non-trivial structure in parameter space, with several secondary minima, sad- 
dle points, valleys, etc. To help the user perform the minimization effectively, 
we have added tools to compute the \ 2 contours and profiles. In addition, 
our program's flexibility allows for many strategies in performing the fit. It 
is therefore possible, by following the techniques described in Sect. 3.7, to 
scan the parameter space and ensure that the minimum found is the true 
one. Further systematic deviations between the model and experiment can 

be recognized via the program's output, which includes a particle-by-particle 

comparison between experiment and theory. 

Purpose: 

In consideration of the wide stream of new data coming out from RHIC, there 
is an on-going activity, with several groups performing analysis of particle 
yields. It is our hope that SHARE will allow to create an analysis standard 
within the community. It can be useful in analyzing the experimental data, 
verifying simple physical assumptions, evaluating expected yields, as well as 
allowing to compare various similar models and programs which are currently 
being used. 

Computation time survey: 

We encounter, in the Fortran version computation, times up to seconds for 
evaluation of particle yields. These rise by up to a factor of 300 in the 
process of minimization and a further factor of a few when xV^doF profiles 
and contours with chemical non-equilibrium are requested. 

Accessibility: 

http : //www. if j . edu.pl/Dept4/ share .html 

or http : // www . physics . arizona/~ torrieri/SHARE/ share . html 
as well as from the authors upon request 

1 Introduction 

In strong interaction reaction processes particle production is abundant. In 
his seminal 1950 article Fermi proposed a non-perturbative description of 
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particle yields based on the assumption that the accessible phase space will 
be fully saturated [3]. Pomeranchuk extended the model by reconsidering 
the particle freeze-out condition [4]. He argued that the reaction volume 
should expand to the point where particles would decouple from each other, 
given their inelastic reaction cross sections. We call this stage the chemical 
freeze-out. In the ensuing decade another important feature of the strong 
interactions was discovered: the existence of numerous hadronic resonances. 
Hagedorn recognized that the large number of different hadronic states is an 
expression and characterization of their strong interactions. He has shown 
that the increase in the number of resonances with their mass has profound 
implications for the behavior of matter at high temperature [5] : hot hadronic 
matter could undergo a phase transition at Hagedorn temperature, Th ~ 
160 MeV. Hagedorn referred to this as the boiling point of hadronic matter 
beyond which the gas of quarks would prevail as a new form of matter. For 
an introduction into the literature and history of this subject see [6], and for 
recent developments also [7]. 

High energy relativistic heavy nuclei (RHI) collisions allow to create, in 
the laboratory environment, a fireball of matter at extreme density and tem- 
perature. The immediate objective of this experimental program is to identify 
formation of a relatively large volume of deconfined quark-gluon matter, the 
quark-gluon plasma (QGP), and to explore its properties. This form of mat- 
ter has existed in the early Universe just before quarks and gluons evolved 
into the 'normal' hadronic particles which surround us today. Formation of 
a region of space in which quarks can roam freely would confirm that their 
confinement is a property of the structured strong interaction vacuum state. 
For a general introduction to these physics questions see [8]. 

Considering the short available lifetime of the fireball, it is difficult to 
develop probes capable of uniquely distinguishing a reaction involving for- 
mation of the deconfined QGP state from the one involving a cascade of 
reactions between individual confined baryons and mesons. The characteris- 
tic feature of the RHI reaction is the formation of a large number of hadronic 
particles. Irrespective if, or not, deconfinement has been achieved, this con- 
stitutes the final state of a RHI collision. Theoretical considerations suggest 
that features of the final hadron abundances are sensitive to the question if 
deconfined state has been formed [10]. 

Analysis of the final hadron state can thus lead to important insights 
regarding reaction mechanisms governing the last stage of the QGP evolution. 
We offer here a standardized version of the Statistical HAdronization with 
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REsonances (SHARE) approach allowing to obtain final state hadron particle 
yields. We have at first prepared two platforms, developed in parallel in the 
form of Fortran and Mathematica shareware, and verified their input and 
logical structure by making sure they yield identical results using the same 
input file. However, in physics applications which require optimization of 
parameters, Mathematica turns out to be forbiddingly slow, and Fortran is 
preferred. For some other applications, the versatility of Mathematica is of 
clear advantage, e.g., one can track the feeding from resonances, make use 
of the integrated graphics, or modify the notebook to extend the model. 

Several research groups beginning with P. Koch et al. [11] have developed 
statistical models of hadron production based on the seminal work by Fermi, 
Pomeranchuk and Hagedorn [3-5]. Hadron yields are computed assuming 
that hadron production can be obtained evaluating the accessible phase space 
size. Known hadron resonances are incorporated and their decay chain eval- 
uated. One has to view this type of hadronization as providing a bottom 
line, there could be always additional microscopic production mechanisms 
which should be most clearly visible for most rarely produced particles. The 
experimental data coming from lower energy collisions of 2-11 GeV at AGS, 
through the CERN SPS energies of 8-17 GeV , to the highest presently avail- 
able energies at RHIC, reaching 200 GeV, can be understood, sometimes to 
a surprising accuracy, with the help of such straightforward statistical ideas. 
There are many variants of statistical approaches, with the common feature 
that at some point of its evolution the hadronic system freezes, and the parti- 
cle fill out the phase space according to statistical distribution. The particles 
form a 'soup' of both stable hadrons, as well as hadronic resonances, which 
later decay, increasing significantly the yields of stable particles. In addition, 
the system expands, i.e., undergoes longitudinal and transverse flow, which is 
an important phenomenon in heavy-ion collisions, distinguishing them from 
elementary hadronic collisions. 

We note that statistical hadronization programs [12-18] require a very 
detailed input of the hadronic spectrum, and definitions of the subsequent 
cascading decays of hadron resonances. Tacit assumptions can make a differ- 
ence of physical significance in the outcome of the analysis. In addition, some 
information is not available for relevant resonances and has to be assessed by 
using the current hadron structure knowledge, in particular regarding parti- 
cle degeneracy and decay patterns. We have excluded from consideration all 
so called single-* resonances, and practically all double-** resonances seen 
in particle data book [2]. The double-** resonances we accept are those 
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where the discovery record is convincing, but a confirmation experiment not 
available. An example are f2-resonances, see comment on f2(2470) in [2]. 

The statistical hadronization approach to describe particle yields arises 
naturally when considering a pot of boiling quark-gluon soup: hadrons evap- 
orate with an abundance corresponding to the accessible phase space. The 
quark-chemical equilibrium in the pot implies that the evaporated hadrons 
are near, but in general not at their chemical equilibrium. For the dynami- 
cally evolving fireball of quark-gluon plasma, a similar situation arises should 
the QGP fireball undergo sudden hadronization in the final stages of its evo- 
lution [19]. In this case QGP breaks rapidly into ever smaller drops of mat- 
ter, ultimately consisting of individual elementary hadronic particle species. 
The statistical particle production model also applies in the opposite limit 
of a very slow hadronization process [20,21], assuming that there is time 
for diverse hadronic particles to be produced and destroyed in a chemical 
(re)equilibration process. 

2 Statistical models in a capsule 

The statistical models we are concerned with have the following main physical 
ingredients: 

• Particle abundances at chemical freeze-out, 

• Resonance decays. 

For the details of different statistical approaches the reader is referred to 
Refs. [22-37]. 



2.1 Statistical particle distributions 

The densities of particle species % are given by the Fermi-Dirac or Bose- 
Einstein distribution functions 

n( m., ffi ;T,T,). n , = 9>J^ r _ leMV J +mim±l , (3) 

- I^M-Tr^A^). (4) 

n=l 

The second form, Eq. (4), expresses the momentum integrals in terms of 
the modified Bessel function K 2 . This form is practical in the numerical 
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calculations and is used in the Fortran code. The series expansion (sum 
over n) converges when Tie~ m ^ T < 1. Consideration of this condition is 
required only for the pion case in the range of parameters of interest. In 
Eq. (4), the upper signs refer to fermions and the lower signs to bosons, 
respectively. Tj is the fugacity factor, and m 8 is the particle mass. The 
quantity = (2Jj + 1) is the spin degeneracy factor as we distinguish all 
particles according to their electrical charge and mass. The index % labels 
different particle species, including hadrons which are stable under strong 
interactions (such as pions, kaons, nucleons or hyperons) and hadron which 
are unstable (p mesons, A(1232), etc.). 

In the most general chemical condition 1 , the fugacity is defined through 
the parameters Xji, X q , A s , A c (expressing, respectively, the isospin, light, strange 
and charm quark fugacity factors), and 7g , 7s , 7 C (expressing the light, strange 
and charm quark phase space occupancies, = 1 for absolute yield equilib- 
rium). The fugacity Tj is then given by: 



T, = A- (A ?7 X* (A.7.)" 2 (Ae7e)^ (A 5 7 5 )^ (A^ (A^f 1 , (5) 



where 
and 



A<? — A ? 1 , A s — A 5 1 , A c — A g 1 , (6) 

lq = lq, Is = Is, 7c = 7c- (7) 

Here, iV*, N l s and N l c are the numbers of light (u,d), strange (s) and charm 
(c) quarks in the ith hadron, and iVi, 7V| and are the numbers of the 
corresponding antiquarks in the same hadron. 

In the case of the models of Refs. [12, 13,29], where the chemical poten- 
tials for the conserved quantum numbers are considered assuming chemical 
equilibrium ( 7? = 7s = 1) and absence of charm (N* — — 0) one has 

T - = exp (BVB + Sys + Ii^ (g) 

where Bi, Si, and Jg are the baryon number, strangeness, and the third 
component of the isospin of the ith particle, and /i's are the correspond- 
ing chemical potentials. In this case, the two formulations are related by 

lr rhis condition is commonly called chemical non-equilibrium. However, the conven- 
tional equilibrium in which existent particles are redistributed according to chemical poten- 
tials is maintained here. The non-equilibrium regarding particle production is, in precise 
terms, called absolute chemical (non)cquilibrium. 
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equations: 

^eq _ e (-3/i s +/i S )/3T^ 

A e 7 ? = X 4 = eW^. (11) 

The user is free to input either format for all of the three chemical potentials 
(fJ>B, Us, f-i>i 3 ) or fugacities (X q , X s , A/ 3 ) (see section 3.1 for details). All quan- 
tities are automatically converted into fugacities for calculations and fits, 
however both chemical potentials and fugacities are presented in the output 
of the program. 

2.2 Resonance decays 

In the first instance, we consider hadronic resonances as if they were particles 
with a given well defined mass, e.g., their decay width is insignificant. All 
hadronic resonances decay rapidly after freeze-out, feeding the stable particle 
abundances. Moreover, heavy resonances may decay in cascades, which are 
implemented in the algorithm where all decays proceed sequentially from the 
heaviest to lightest particles. As a consequence, the light particles obtain 
contributions from the heavier particles, which have the form 

ni = b 2 ^\...b N ^ rN - 1 n N , (12) 

where bk^k-i combines the branching ratio for the k — > k—1 decay (appearing 
in [2]) with the appropriate Clebsch-Gordan coefficient. The latter accounts 
for the isospin symmetry in strong decays and allows us to treat separately 
different charged states of isospin multiplets of particles such as nucleons, 
Deltas, pions, kaons, etc. For example, different isospin multiplet member 
states of A decay according to the following pattern: 

A ++ ^7T + +p, (13) 

A+^i(7r+ + n) + ^(7r°+p), (14) 

A°^^-+p) + ^ + n), (15) 
A~^7r~+n. (16) 
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Here, the branching ratio is 1 but the Clebsch-Gordan coefficients introduce 
another factor leading to the effective branching ratios of 1/3 or 2/3, where 
appropriate. 

For two-body strong decays (such as the example shown above), the 
Clebsch-Gordan coefficients are automatically calculated by the program. 
In this case, the branching ratio taken from [2] is the only input in the 
program. For three-body decays, the products of the Clebsch-Gordan coef- 
ficients have been averaged, with equal weights, over all possible couplings 
between the three isospins. The result of this procedure has been used to 
rescale the branching ratio taken from [2], and the effective branching ratio 
has been obtained in this way, which is delivered as the input. 

The case of A — > ttN is easy to deal with, since only one decay channel 
is present and the branching ratio is well known. However, in most of the 
cases, we take into account several decay channels appear. We note that 
the partial widths (product of branching ratio with total width) are not 
sufficiently well known. In our approach we disregard, as a rule, all decays 
with the branching ratios smaller than 1%. In addition, if the decay channels 
are classified as dominant, large, seen, or possibly seen, we always take into 
account the most important channel. If two or more channels are described 
as equally important, we take all of them with the same weight. For example 
/o(980) decays into nil (according to [2] this is the dominant channel) and 
KK (according to [2] this is the seen channel). In our approach, according to 
the rules stated above, we include only the process /o(980) — ► Tin. Similarly, 
ao(1450) has three decay channels: r/n (seen), 7r?/(958) (seen), and KK 
(again seen). In this case, we include all three decay channels with the 
weight (branching ratio) 1/3. Certainly this procedure is not unique and 
other methods of selection/treatment of the decays are conceivable. 

Another difficulty is that the branching ratios are not given exactly (in- 
stead of one value, the whole range of acceptable values is given) and the sum 
of the branching ratios may differ significantly from 1. In this case we take 
the mean values of the branching ratios. Since we require that their sum is 
properly normalized, sometimes we are forced to rescale all the mean values 
in such a way that their sum is indeed 1. Of course, a different convention 
from ours can be implemented by the user through a different particle data 
input file. 

Weak decays (identified by the fact that they break strangeness and 
isospin) are stored separately, and added to the ratios with a feed-down cor- 
rection set by the user. However, electromagnetic decays such as S° — > A + 7 
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are counted as if it were hadronic in the contribution to the yield of A (and 
the antiparticle). 

2.3 Finite resonance width 

If the particle % has a finite width Tj, the thermal yield of the particle is more 
appropriately obtained by weighting Eq. (3) over a range of masses to take 
the mass spread into account: 



The use of the Breit-Wigner distribution with energy independent width 
means that there is a finite probability that the resonance would be formed at 
unrealistically small masses. Since the weight involves a thermal distribution 
n(M, gf, T, Tj) which would contribute in this unphysical domain, one has to 
use, in Eq. (17), an energy dependent width. 

The dominant energy dependence of the width is due to the decay thresh- 
old energy phase space factor, dependent on the angular momentum present 
in the decay. The explicit form can be seen in the corresponding reverse 
production cross sections [38,39]. The energy dependent partial width in the 
channel i — > j is to a good approximation: 



Here, is the threshold of the decay reaction with branching ratio b^j. For 
example for the decay of i :— A ++ into j := p + 7i + , we have = m p + m n +, 
while the branching ratio is unity and the angular momentum released in 
decay is = 1. From these partial widths the total energy dependent width 
arises, 






(18) 




(19) 



j 



For a resonance with width, we thus have replacing Eq. (17): 



n, 



^ = ^E / dMn(M, gi] T,T t ) 



(20) 



(M-m^ + tr^M)] 2 /^ 
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and the factor N (replacing 2n) ensures the normalization: 

N- = V f°dM r ^' (M) (21) 

In principle, Eq. (20) does not take into account the possibility that the 
state into which one is decaying is itself a unstable state in a thermal bath. 
Doing this would require a further average over the width distribution of 
receiving state. We will not implement here this refinement. 

The evaluation of the Breit-Wigner integral requires additional attention. 
The presence of the thermal distribution limits the integral range in Eq. (17) 
and thus the integral is finite and well behaved. However, it involves in gen- 
eral diverse cases of both very narrow, r<m, and wide, T ~ m, resonances. 
To perform the integral reliably and precisely in all of these cases we split it 
in two parts, = f™° +2r + J™ +2V ^ ne ^ rs ^ P ar ^ can ^ e securer y done 
using, e.g., Gaussian integration. The second part, with a potentially slowly 
falling tail, can be computed through a variable change z = < - m ° A ^ 2r ' ) before 
the numerical integration. 

Calculating widths adds considerably to the computation time. It should 
be noted that Eq. 3 is in the form of K 2 functions depending on temperature 
and mass, weighted by coefficients representing chemical potentials. The 
expensive width integrals, therefore, only have to be recalculated when the 
temperature under consideration changes. The program takes advantage of 
this in its fitting algorithm, something which can dramatically increase fitting 
time. See section 3.7 for details. 



3 SHARE structure 

SHARE is a modular program, whose basic structure is illustrated in Fig. 1. 
Five input data files are required, containing model and experimental data. 
There is a calculational and a fitting block, piloted by instructions which the 
program reads from a running file (called sharerun.dat a). Each command 
will generally be performed independently from the others, and generate its 
own output in a separate file. 

Any user defined file names of specified length can be used for all input 
files except sharerun.data, but the defaults are: 

thermo.data (11 letter filename), 



12 



Thermal 

parameters 

thermo.data 



Particles 
Data 

particles. data 



Resonances 
decay tree 
decays.data 



Statistical hadronization 
sharerun.data 

calc totratios 
calc ratioplot 
calc ratiocont 



Particle ratios 
Extensive 
quantities 
Calculation 
output file 
(13-letters) 



Experimental 
data 

totratios. data 





Fitting routines 



calc fitratios 
calc chiprofile 
calc chi2_cont 



Best fit parameters 

parameter errors 

statistical significance 

point by point model-data comparison 

Fitting output file (13 letters) 

parameter sensitivity 

number of minima 

correlations 

% profiles, contours file (12 letters) 



Figure 1: Structure of the SHARE package. Running commands which can 
be given are in red, default input filenames are in blue, while possible output 
files are in violet (color online). All output filenames are set by the user 
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particles . data (14 letter filename), 
decays. data (11 letter filename), 
totratios . data (14 letter filename), 
ratioset .data (13 letter filename). 

These deal, respectively, with initial thermal parameters defining the be- 
ginning of the fit procedure, particle properties, decay patterns, and the par- 
ticle ratios (Math SHARE uses particles . data and decays, data only). 
These input files are described in detail in the next subsections. It is pos- 
sible to insert comments into all of these files: any time the first character 
starts with the subsequent input is disregarded (CAUTION: For Math 
SHARE a comment line must include at least one more word after e.g., 

This is a comment' ). 

3.1 Thermal parameter input 

(11 letter filename, default name: thermo.data) 

The thermo . data file contains input temperature and chemical potentials. 
Where applicable, all units are GeV. A sample is given in Table 1. 



tag 


value 


explanation (not part of file) 


temp 


0.165 


temperature 


mu_b 


0.028 


light quark fugacity or chemical potential 


muj3 


0.006 


strange quark fugacity or chemical potential 


gamq 


1. 


light quark phase space occupancy 


gams 


1. 


strange quark phase space occupancy 


lmi3 


-0.001 


Is fugacity or chemical potential 


norm 


1. 


absolute normalization 


lame 


1. 


charm quark fugacity 


game 


0.001 


charm quark phase space occupancy 


accu 


0.001 


calculations accuracy 



Table 1: A typical thermal input file 



Chemical potentials can be input as X's and 7's, Eqs. (5)-(7), or as fjfs, 
Eqs. (8)-(10). If the tag is, respectively, 'mu_b', 'mu_s' or 'mui3' the program 
assumes input is given as in Eqs. (5)-(7) otherwise Eqs. (8)-(10) are assumed. 
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The 4-letter tags given in the file are used throughout the program to label 
the corresponding thermodynamic parameters (fitting variables, arguments 
for profiles and contour plots, etc). 

3.2 Particle properties data file 

(14 letter filename, default name: particles. data) 

This input file contains information about the properties of particles such 
as: mass, width, spin, isospin, the quark contents, and the Monte-Carlo 
identification number. The data file is written in the following format 

name mass width spin I 13 q s aq as c ac MC 

where: 

name - a nine-letter character string identifying the particle, 

mass — mass in GeV, 

width — width in GeV, 

spin — spin, 

I — isospin, 

13 — 3rd component of isospin, 

q, s — number of light /strange quarks, 

aq, as - number of light /strange antiquarks, 

c, ac - number of charm/anticharm quarks, 

MC - particle's identification number, usually (where applicable) cor- 
responding to the standard Monte-Carlo particle identification conven- 
tion [2]. 

For example, the A(1232) ++ will appear in the input file as 
D11232plp 1.2320000 0.1200000 1.5 1.5 1.5 3 2224 

The quark number can be a non-integer, to accommodate strong interac- 
tion flavor mixing such as that of the rj. Note that SHARE calculations are 
relevant for a strongly interacting system, where the relevant states are 
and K°. K° — K° mixing is an electroweak process, occurring on a longer 
timescale, and should be implemented at the end of the calculation. 

Several versions of this input file are available, described in detail in section 
6.1 
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Our particle naming convention is to form a 9-letter name through a letters- 
mass-ending combination (e.g., Lmlll5zer for A), with the following usual 
three letter endings 

zer for Zero 
zrb for Zerobar 
plu for plus 
plb for plusbar 
min for minus 
mnb for minusbar 
pip for plusplus 
ppb for plusplusbar 

sht for particle with this ending sht = ^= (zer+zrb) 

e.g., K0492sht for K s = 75(^0 + #0) 
tot tot = (zer+zrb) 

e.g., Lmlll5tot for A + A 
mnt mnt = (min+mnb) 

e.g., UM1321mnt for fl +H 
pmb pmb = (plu+plb) 

e.g., pr0938pmb for p + p 
plm plm = (plu+min) 

e.g., Xil321plm for 2+ + S" 

These conventions will be assumed also when fitting particle ratios and/or 
yields (see section 3.4 for details) 

3.3 Particle decay pattern data file 

(11 letter filename, default name: decays. data) 

This input file contains the information on particle decays. The data file 
lines are here written in the format: 

Name parcnt Name daughtcrl Name daughtcr2 Name daughtcr3 2 BR C-G?(0/1) 

where BR refers to the branching ratio of this decay (without the Clebsch- 
Gordan factor, as it appears in the data book) and C-G refers to whether 

2 Where applicable. 
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the branching ratio should be completed by a Clebsch-Gordan coefficients 
(0: no, 1: yes). For instance, the decays A + — > 7r + + n and A + — > 7r° +p will 
be given as (compare Eqs. (13)— (16)): 

D11232plu pi0139plu ne0939zer 1.0 1 
D11232plu pi0135zer pr0938plu 1.0 1 

while f] decays will be (see discussion in Sect. 2.2): 

eta547zer gamOOOzer gamOOOzer 0.3943 

eta547zer pi0135zer pi0135zer pi0135zer 0.3251 

eta547zer pi0139plu pi0139min pi0135zer 0.226 

eta547zer pi0139plu pi0139min gamOOOzer 0.0468 

3.4 (Experimental) Values to be calculated 

(14-letter filename, default file: totratios.data) 

Experimental data values of interest are submitted in in the following format. 



3.4.1 Particle ratios 

When considering the particle ratio evaluation: namel name2 data ran- 
dom systematic fit? (-1/0/1/2) 

where 

namel The first particle in the ratio (numerator) 

name2 The second particle in the ratio (denominator. Can also be a 

tag indicating the quantity is not a ratio but a yield or a density) 
data The experimental value of the data point 

random The random (statistical) error 
systematic The systematic error 

fit? This ratio contributes to the evaluation of x 2 /Nd f if this pa- 

rameter is set to 1 or 2. If the parameter is set to 0, the ratio is not 
fitted, but calculated and output to the graph file (see cards in section 
4 for details), if fit = —1 or 2 means the ratio is not output to the 
graph file. 



17 



Strangeness conservation 



# 

# 

# Exact 



1r 

netstrang 


solvelams 


0. 


0.0001 


0. 




2 


# Approximate 












#netstrang 
# 


totstrang 


0. 


0.0001 


0. 




2 


#Weak decay scheme 


K_S->pi 


K_L->pi 


Y->pi, K 


Y- 


->N 


weakdecay 


top 


0.0 


0.00 


0.0 


1 


. 00 


weakdecay 
# 


bottom 


0.0 


0.00 


0.0 


1 


.00 


#ratio 


ratio 












♦numerator 
# 


denominator 


yield 


star, err 


sys . err. 




Fit? 


Xil321min 


pi0139min 


0.0093 


0.002 


0. 




1 


Xil321mnb 


pi0139min 


0.0080 


0.0017 


0. 




1 


Xil321mnb 


Xil321min 


0.8530 


0.1000 


0. 




1 


UM1672mnt 


Xil321mnt 


0.1438 


0.035 


0. 




1 


Xil321min 


Lmlll5zer 


0.1870 


0.460 


0. 




1 


Xil321mnb 


Lmlll5zrb 


0.2150 


0.5400 


0. 




1 


Xil321min 


phl020zer 


0.3400 


0.070 


0. 




1 


Lmlll5zrb 


Lmlll5zer 


. 74 


0.05 


0. 




1 


Ka0492min 


pi0139min 


0.181 


0.025 


0. 




1 


Ka0492plu 
# 


pi0139plu 


0.201 


0.025 


0. 




1 


#Weak decay scheme 


K_S->pi 


K_L->pi 


Y->pi, K 


Y- 


->N 


weakdecay 


top 


0.5 


0.00 


0.5 


1 


.00 


weakdecay 
# 


bottom 


0.5 


0.00 


0.5 


1 


.00 



Ka0492plu 
Ka0492plu 

# 

pr0938plb 

# 

#STAR 
Ka0492plu 
#PHENIX 
Ka0492plu 

# 

# 

Lml520zer 
rho770zer 
f 00980zer 
Ka0892zer 



pr0938plu 
pr0938plb 



0. 
0. 



pr0938plu 



Ka0492min 
Ka0492min 



2.4 0.5 

3.0 0.5 

"other ratios" 

0.71 0.06 0. 
K+/K- Both experiments 

1.103 0.025 0. 

1.153 0.040 0. 



Lmlll5zer 
pi0139min 
pi0139min 
Ka0492min 



RESONANCES 
. 022 
0.1 
0.02 
0.24 



. 01 
0.02 
0.01 
0.01 



0. 
0. 
0. 
0. 



#Weak decay scheme 


K„S->pi K_ 


_L->pi 


Y->pi, K 


Y->N 


weakdecay 


top 


0.50 


0.00 


1.00 


0.00 


weakdecay 
# 


bottom 


0.50 


0.00 


1.00 


0.00 


phl020zer 


Ka0892min 


.595 


0.24 


0. 


1 


D11232plp 


pr0938plu 


0.24 


0.037 


0. 


1 


# Predictions (NOT fitted) 


#Sgl385plu 


Xil321min 


. 64403 


0.01 


0. 





#Xil530min 


Xil321min 


. 1802 


. 01 


0. 





tpentaqOth 


Xil530min 


1 . 1437 


0.10 


0. 





#rho770zer 


phl020zer 


. 1 


0.02 


0. 





#f 00980zer 


phl020zer 


0.1 


0.02 


0. 






Figure 2: A typical totratios.data file 



3.4.2 Particle yields and source properties 

The entries in this files have a second use when namel and name2 are not 
particle names but are as follows: 
In case name2 is: 

prt_yield the yield of the first particle or collective quantity. 

prdensity the density of the first particle or collective quantity (in fm~ 3 ). 

solveXXXX See section 3.6 

In case namel is: 

negatives All negative particles stable under the strong interaction 
totstrang strangeness 3 (s) 
netstrang net strangeness (s — s) 
totenergy energy (in GeV) 

totbaryon sum of all baryons and antibaryons, B + B 

netbaryon net baryon number, i.e., baryons minus antibaryons, B — B 

totcharge charge Q 

net charge net charge < Q — Q > 

entropy _t entropy S 

We provide several data files as an example, described in section 6.1 

3.4.3 Weak decay feed-downs 

All weak decays are included with the provided SHARE decay tree input 
files. They are identified as weak by the program due to the violation of 
strangeness and isospin within the decay. Weak decays present additional 
complications, since their reconstruction acceptance is usually non-trivial 
and experiment-specific. 

SHARE therefore, includes additional parameters, set by the user, to 
regulate the acceptance of these decays. These parameters are read in from 
the same file as the experimental data points, and have approximately the 
same format. 



3 For this and the subsequent quantities, volume normalization can be correct only when 
total particle yields are considered. The units of the calculated result vary according to 
whether the quantity is a density, a yield or a ratio 
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A statement specifying weak decay feed-downs has the form 
weakdecay Tag F 1 F 2 F 3 F 4 
Where: 

Tag specifies weather the coefficients refer to the numerator or denominator 
of the ratio (note that they are often different, e.g. in a ratio such as 
S//i~). If the user is setting the numerator feed-down coefficients, tag 
is set to "top" . In case of the denominator, tag is set to "bottom" . 

Fi_ 4 are four numbers, specifying the weak feed-down coefficients. 

• Fi refers to the acceptance of re from K s — > decays. It is by 
default set to 0.5 

• F 2 refers to the acceptance of n from K L — > nnn decays. It is 
typically very small, so it is by default set to 

• F 3 refers to the acceptance of mesons from weak hyperon decays. 
It is by default set to 0.25 

• F4 refers to the acceptance of baryons from weak hyperon decays. 
It is by default set to 1. 

It is possible to specify a different feed-down for each data-point, by preceding 
a data line with "weakdecay" statements, eg: 

weakdecay top 1.0 0.0 0.3 0.8 

weakdecay bottom 0.0 0.0 0.5 0.8 

Xil321mnt negatives 0.2 0.05 0. 1 
It is also possible to specify one weak feed-down for a series of data points, 
or all the datapoints analyzed; The program assumes all data points after 
a "weakdecay" statement to have the same value until a new "weakdecay" 
statement is encountered. For example, 



weakdecay 


top 


1.0 


0.0 


0.3 


0.8 


weakdecay 


bottom 


0.0 


0.0 


0.5 


0.8 


Xil321mnt 


negatives 


0.2 


0.05 


0. 


1 


Lmlll5zer 


pr0938plu 


0.15 


0.05 


0.1 


1 



fixes the feed-down for both E/h~ and A/p ratios. (In this case, of course, 
three of the four feed-down coefficients will have no effect on the second 
data-point) 
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3.5 Fit parameters 

(13 letter filename, default name: ratioset.data) 



The ratioset.data file defines the parameters which will be varied during 
the fit or kept constant. It further contains the limits to be used in the fit 
(for the use of MINUIT), as well as the initial step size. The typical format 
is: 

tag lower / upper limit step size fit? (0/1) 



temp 


0. 


1. 


0.01 


1 


lamq 


0. 


10. 


0.1 


1 


lams 


0. 


10. 


0.1 


1 


gamq 


0. 


10. 


0.01 


1 


gams 


0. 


10. 


0.1 


1 


mui3 


0. 


10. 


0.2 


1 


norm 


0. 


10000. 


0.3 





lame 


0. 


10. 


0.1 





game 


0. 


10. 


0.1 






If the lower and upper limit are equal, MINUIT will fit with no limits. If 
the only fit parameters are particle ratios and densities, the normalization 
is automatically kept fixed. The parameters can be input in any order. We 
have found that the fit quality (speed, reliability) depends considerably on 
the order of parameter input which is retained calling the MINUIT package: 
firstly, the most significant fit parameters should be input first (temp, norm). 
Secondly, the highly correlated parameters should be placed next to each 
other. 

Note that values of parameters also arise as result of conservation laws 
rather than from fits to particle yields. For example, A s can be computed 
in terms of the other parameters requiring overall strangeness conservation 
within the phase space domain covered by the results considered. A/ 3 is sim- 
ilarly determined via the participant matter proton-neutron ratio, however, 
here it is possible to use this only when the data is available in 4n accep- 
tance, since other hadrons can buffer the balance of charge condition. For 
example, at RHIC in the central rapidity region the large number of pions 
implies in effect nearly symmetric A/ 3 = 1. Also, /i# can be fixed through 
the participant number in the phase space region observed. The next section 
describes how SHARE allows the user to implement these constraints either 
by solving the constraints (numerically) or by a fit. 
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3.6 Implementing conservation laws 

SHARE allows the user to solve for a fit parameter, rather than to find its 
value through a fit. For instance, it is possible to implement strangeness 
conservation by solving numerically the constraint equation 

(s-s)=0 (22) 

for A s . In this case, A s is not a fit parameter anymore, but rather an analytical 
function of the other fit parameters and experimental particle yields. 

To solve for a fit parameter, name2 in the totratios.data file should be of 
the form solveXXXX where XXXX corresponds to the parameter for which 
one wants to solve. The parameter limits are still used by the program, as a 
constraint solution outside the limits is rejected. This is useful for rejecting 
unphysical solutions of the constraint equation, such as X s < 0. 

It is, in principle, possible to solve any data point for any thermal param- 
eter. However, many such combinations do not have minima to which the 
equations converge nicely. If this is the case with a lot of MINUIT iterations, 
it is unlikely that the minimization procedure will work. 

It is therefore recommended that the solving algorithm only be used to 
solve for chemical potentials from conservation laws. For instance, to ensure 
strangeness conservation through solving for \ s the input file should contain 
the following line: 

netstrang solvelams 0. 0. 0. 0. 0. 

A line solving for X q from baryon conservation in case of Pb-Pb collisions 
might be 

netbaryon solvelamq 362. 0. 0. 0. 0. 

and the corresponding charge conservation statement will be 
netcharge solvelmi3 142. 0. 0. 0. 0. 

Note that it is also possible to implement conservation laws in terms of 

particle ratios. The two lines 

netbaryon solvenorm 362. 0. 0. 0. 0. 

netcharge solvelamq 142. 0. 0. 0. 0. 

will fix \ q in terms of the baryon/charge = ^p- = 0.44 for Pb-Pb) ratio 
even when the absolute normalization norm is a dummy variable which 
appears in no other data point and is not used in the fit. 

The alternative to exact solving is to implement a conservation law by 
treating it as a data point. A line such as 
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netbaryon prt_yield 362. 10. 0. 0. 0. 1 

will make sure that the baryon number is close to the one expected for Pb-Pb 

collisions at SPS. Similarly, 

netstrang totstrang 0. 0.01 0. 0. 0. 1 

will make conserve strangeness to one unit in 100 s pairs, rather than solving 
the constraint exactly. 

The choice of whether to implement the conservation law analytically or 
through a fit is therefore left to the user. It was found that both approaches 
are giving very similar results. However the most reliable approach is to use 
the exact constraint only. 

solve statements should be put at the top of the experimental ratios file. 
If this is not done, the program returns with an error. 

3.7 Fitting strategies 

Multi-parameter fits can be very time-consuming, and occasionally give the 
wrong answer when a "false" (i.e., local rather than absolute) minimum is 
found. If these minima are close enough, MINUIT's error and contour cal- 
culations fail to converge convincingly. If this is the case, a warning message 
is printed in the minimization output file. 

The user should make note of a few simple rules which can dramatically 
decrease the fitting time, and enhance the probability that the right minimum 
is found. 

• Appropriate calibration of the step size and limits can: 

- dramatically decrease the running time for the case that the user 
already knows a good estimate of the minimum, 

- select the appropriate minimum if several minimal points are present 
in the parameter space. Certain datasets analyzed in the litera- 
ture present three minimal points: at equilibrium, at 7 g > 1, and 
lq < 1 [18]. 

• Fitting absolute yields is more difficult (prone to converge to a false 
minimal points or spend very long time finding a minimum), as the 
yield normalization strongly correlates with other parameters of the 
fit. Especially with a limited experimental data set it is not always 
easy to distinguish a hot, dense and small fireball from a colder, more 
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dilute larger one using yields alone. It is therefore advisable to convert 
all data points to ratios. 

• Inclusion in the fit of extensive fireball properties, which require a sum 
over every hadron, considerably increase the running time. This is es- 
pecially true for energy (density), pressure and entropy (density), which 
require the Bessel function evaluation to be redone for every hadron at 
each iteration step. It is therefore suggested that the experimental file 
used for the fit not contain these quantities: They can be calculated 
after the fit in a single pass evaluation, with the fitted parameters as 
input, using modified input files. The provided experimental input files 
(see section 6.1) contain the thermodynamic quantities used in a typi- 
cal fit (total baryon number, strangeness and charge). We also provide 
an input file (totratext.data) devoted exclusively to thermodynamic 
quantities. 

• The more variable parameters a fit has, the more minimal points are 
likely. In case of a particularly difficult fit, it is suggested it be com- 
pleted in two stages: An initial fit with only temperature and normal- 
ization as varying parameters, followed by a subsequent fit, starting 
from the previously found minimum, where the remaining parameters 
are released (the sharerun.data provided with the program is based 
on this strategy). Alternatively, the conservation law constraints (for 
e.g., zero net strangeness, given charge and/or baryon number) can 
be tightened or relaxed. A x 2 profile can also be helpful. Generally, 
a minimum with a x 2 significantly above expectation signals that the 
true minimum has not been found. 

• Fits with finite resonance widths are much more computationally ex- 
pensive. However, as explained in section 2.3, the computationally 
expensive integrals and Bessel functions only have to be recalculated 
when the temperature is changed. For this reason, a temperature x 2 
profile (which finds the minimum at each value of the temperature pa- 
rameter, kept fixed throughout each fit. Hence, the temperature only 
has to be calculated once per datapoint) can be computed much more 
quickly than a fit in which temperature is a free parameter. 
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4 Running program cards (sharerun.data) 



This file contains the instructions which the program executes on running. 
Each line corresponds to a different operation, such as reading data files, 
assigning values to parameters, calculating ratios, minimizing, and plotting 
contours and x 2 /NDoF profiles. The program will read this file line by line, 
execute each command, and stop when it reaches the end of the file. Fig. 1 
has an example of a typical sharerun.data file. 4 

Each command can be used more than once with different input and 
output files, one at a time. We shall proceed to give a detailed description 
of each command's meaning and syntax. Keep in mind that 2 spaces have 
to be maintained between each word or number. The lines have the general 
form: 

PSET <4-letter tag> VALUE 

This Parameter SET command sets the thermodynamic variable de- 
fined by the TAG (name) to its designated VALUE. The initial com- 
mand shown in figure 3, READ THERMJNI, comprises a series of 
PSET-type commands which read from an input file, covering all al- 
lowed thermal parameters. 

READ THERMJNI <ll-letter filename> 

Reads the file corresponding to thermo.data, containing the values for 
the thermal parameters, as described in section 3.1. 

READ PARTICLES <14-letter filename> <ll-letter filename> 

Reads the file containing particle properties as well as the file containing 
the decay tree, as described in sections 3.2 and 3.3. 

READ TOTALDATA <14-letter filename> 

Reads the file containing experimental data and particle ratios to cal- 
culate, as described in section 3.4. 

PFIX tag Fixes the given thermodynamic Parameter to its current value. 
In a fit, that parameter will not be a variable. 

PREL tag <Lower limit> <Upper limit> <Step size> Releases the 
given thermodynamic Parameter, giving the limits and initial step size 
used in the fit. The next command, 'READ FIT_PARAM', can be 

4 Thc program, being written in FORTRAN, is very sensitive to the file's format. An 
extra space may cause the sharerun.data file to be unreadable. Hence, it is recommended 
that the user bases his modifications on the sharerun.data file shown in the Fig. 1. 
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# Reading experimental data 
READ TOTALDATA totratios . data 

# no width (file partnowdt . data) , should not take that long 

READ PARTICLES partnowdt . data decays. data 

# Temperature initialization file (CAN ALSO USE PSET, as in width calculation) 



READ 


THERM_INI 


thermo . data 








READ 


F I T_P ARAM 


rtonlnrm . data 








CALC 


FITRATIOS 


f it_ junkl 


. out 








READ 


F I T_P ARAM 


ratioset . data 








CALC 


FITRATIOS 


f it_nwrhi 


. out 








CALC 


PLOT_DATA 


graph_nwd 


.fit graph_ 


_nwd . 


clc 


graph_nwd 


READ 


TOTALDATA 


totratext 


. data 








CALC 


RAT I OD ATA 


extrhicnw 


. out 








READ 


TOTALDATA 


totratios 


. data 








CALC 


SIGPROr 1L 


prof_T_nowdt temp 


. 12 


. 1 


8 100 


CALC 


SIGPROFIL 


profgq_nowdt gamq 


. 12 


1 . 8 


50 


PSET 


gamq 1 . 












PSET 


gams 1 . 












PFIX 


gamq 












PFIX 


gams 












CALC 


FITRATIOS 


f it_nweqi 


. out 








CALC 


PLOT_DATA 


graph_neq 


.fit graph_ 


_neq. 


clc 


graph_neq 


READ 


TOTALDATA 


totratext 


. data 








CALC 


RAT I OD ATA 


extrhicne 


. out 








READ 


TOTALDATA 


totratios 


. data 








CALC 


SIGPROFIL 


prof_T_noweq temp 


. 12 


. 1 


8 100 


PREL 


gams . 


1 5. 


0.02 








CALC 


FITRATIOS 


f it_rhins 


. out 








CALC 


PLOT_DATA 


grphrhins 


.fit grphrhins. 


clc 


grphrhins 


READ 


TOTALDATA 


totratext 


. data 








CALC 


RAT I OD ATA 


ext_rhins 


. out 








READ 


TOTALDATA 


totratios 


. data 








CALC 


SIGPROFIL 


prof_T_rhins temp 


0.09 


. 1 


8 100 


CALC 


SIGPROFIL 


prof gs„rhins gams 


0.03 


1 . 8 


50 
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Figure 3: A typical sharerun.data file 



understood as a series of PFIX and PREL instructions read from an 
input file. 

READ FIT_PARAM < 14- letter filename> 

Reads the file containing fit parameters, as described in section 3.5. 
CALC RATIODATA <13-letter filename> 

Calculates the value of the ratios read with the READ RATIODATA 
command and the current value of the thermal parameters (obtained 
either from READ THERM_INI or a previous fit). The output of the 
calculation is stored with the given filename, as a table in the following 
format: 

RATIO NAM1/NAM2 <numerical value> 
CALC RATIOPLOT Datapoint <12-letter Filename> tag L H P 

Calculates the ratio of 2 particles, or a particle's yield, or a thermody- 
namic quantity, as a function of the variable represented by the 4-letter 
tag, in a parameter segment delimited by the limits (L,H) and number 
of points (P). The ratio to be calculated must be in the experimental 
ratios datafile; datapoint specifies which point to graph (for instance, 
if datapoint is 2, the second point from the to will be calculated). The 
output is saved as a 2-column table in the given output file, and can 
be plotted with packages such as PAW, Mongo or Xmgrace. 
CALC RATIOCONT datapoint <12-letter Filename> 
tagl LI HI PI tag2 L2 H2 P2 5 

Calculates the ratio of 2 particles, or a particle's yield or density, or a 
thermodynamic quantity, as a function of two thermodynamic variables 
represented by the two 4-letter tags, in the parameter region delimited 
by the limits (LI, 2 HI, 2) and number of points. (PI, 2. Note that an 
100 x 100 grid has 10 000 points, and can take a long time to calculate). 
The ratio to be calculated is indicated in the same way as in CALC 
RATIOPLOT. The output is saved as a 3-column table in the given 
output file, and can be plotted with a program capable of 3D plots. 
CALC FITRATIOS <13-letter filename> 

Minimizes the xV^doF of the set of experimental data obtained through 
READ RATIODATA according to the parameters in 
READ FIT_PARAM and initial values in 

READ THERM JNI. The output is written out to the given filename 
in the following format: 

5 one line in the file 
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First, the output parameters (+/— error if fitted). 

Then the detailed fit results, as a table with the format: 

Top Bottom Theory Experiment Tot. error. Chiterm 

Where Top and Bottom refer to each ratio's numerator and denomina- 
tor, and Chiterm for each ratio refers to x as defined in Eq. 1 

x = /experiment -/theory (= Q . f ^ ^ 
^/statistical ^/systematic 

Finally, the total x 2 /N DoF is presented. A typical output file is shown 
in Fig. 2 

CALC PLOT_DATA <3 13-letter filenames> 

Generates three files which are optimized to be graphed by a package 
such as PAW, Mongo or Xmgrace. The first file has a numerical list of 
ratios which were fitted, the second, a numerical list of calculated, but 
not fitted, ratios. The third has the experimental data, including the 
error bars. See the discussion at the bottom of section 3.4 for details 
on which ratio gets fitted and which just gets calculated. 

CALC CHIPROFIL <12-letter file> tag L H P 

This computes a xV^doF profile of the Parameter designated by tag 
<Parameter tag> (see the thermal input file for a list), from the L to 
the H limit (real numbers), with P specifying the number of computed 
points (integer). The given file will store the main result, as a 2 column 
table of the parameter value and xV^doF- The minimum of each of 
the other parameters for each data point will also be stored in files in 
which the parameter is appended to the name. For instance, if the T, 
X 2 /NdoF profile is stored in file 'profTforratio' the minimal points of 7 g 
across the xV^doF profile are stored in file 'profTforratio_gamq'. 

CALC SIGPROFIL This command is very similar to 

CALC CHIPROFIL. However, instead of the x 2 , the profile for the 
statistical significance is calculated. The command syntax is identical 
to CALC CHIPROFIL 

CALC CHI2_CONT <9-letter filename> deviation tagl tag2 Computes 
a x 2 /NdoF contour of parameters denoted by tagl and tag2, with a 
given deviation from the x 2 /DoF minimum (e.g., 1.1 for a 1.1 x Xmin 
contour). The contour is stored in an output file 9 characters long. 
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4.1 Run log (sharerun.out) 

The complete 'log' for each run is saved in the file sharerun.out. This includes: 

• The content of each input file (in the same format as read) 

• A list of performed operations 

• The output from MINUIT 

• The content of each output file (in the same format as in the file) 

If the program ends without a problem, the message 'PROGRAM TERMI- 
NATED SUCCESSFULLY' is printed, both on the screen and this file. If an 
error occurs, the program writes to the screen that an error has occurred, 
and outputs the details of the error to sharerun.log. 

5 Installation 

SHARE is distributed in a form of an archive containing source and data files. 
Running under the Linux operating system is supported presently. We have 
checked SHARE on FORTRAN in a few Linux installations. In order to run 
SHARE one needs: 

• FORTRAN compiler e.g. g77, f77; 

• CERN library 

SHARE needs the CERN Library, in particular MINUIT and the CERN 
mathematical routines library. These can be found on the CERN web- 
site if not installed in current environment. 

After unpacking, one needs to compile SHARE using the included shellscript 
'fortrat' (which needs to be amended with the correct path to where the 
CERN libraries are) 

6 Organization of the Fortran package 

6.1 Directory tree 
source code: 

sharevl.l.f F77 
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sharel.O.nb mathematica 
Data files: 

thermo.data Thermal inputs 

particles. data Particle datafile — full width 

partnowdt.data Particle datafile — no width 

decays. data Full decays file 

dec_no.data Empty decays file 

ratioset.data Minimization settings 

rtonlnrm.data Minimization settings with most parameters fixed (used in 
the course of an intermediate minimization, see section 3.7) 

totratios.data Sample input drawn from RHIC data (2004). 

totratext.data Sample input with extensive thermodynamic functions 

chi2.data Sample input in Mathematica-readable format 

sharerun.data running script optimized for RHIC 

Shell-scripts: 

fortrat Sample compiling shell-script 
tarrat Sample archiving shell-script 

7 SHARE webpage edition 

A reduced version of SHARE is also accessible directly from the web, at the 
same website from which the program is downloaded. It is possible to submit 
thermal parameters, and calculate a select table of experimentally relevant 
ratios and thermodynamic quantities. The webpage can be used as a cross- 
check between statistical hadronization packages, and for a comparison with 
experimental data. For full input and calculational versatility, however, we 
recommend downloading the full version 
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8 Mathematica version of SHARE 

The Mathematica notebook Math SHARE, included in the SHARE distribu- 
tion, has a similar functionality as the Fortran program but greatly restricted 
physical scope. It uses the same universal input files as the Fortran, namely 
particles . data, decays .data, and in addition the file chi2.data. The 
notebook is designed to perform three types of calculations: 

1. Run the thermal model for specific thermal parameters T, /i B , fi s , and 
fij (in the Mathematica code the 7 factors of Sect. 2.1 are set to unity) 
and determine the particle ratios, 

2. Fit the thermal parameters to the experimental ratios read from file 
chi2.data via the \ 2 method, 

3. Track the feeding of a specified particle from decays of higher states. 
This shows the physical importance of high-lying states and allows to 
determine the amount of feeding for a given state. 

The installation of Math SHARE requires only setting the path to the 
data files (the first cell of the notebook). The following modules are partic- 
ularly useful when running the notebook: 

• readpart [f n] reads the particle properties from the file f n and initial- 
izes the primordial yields to zero 

• decay [fnl] performs the decays according to the file fnl (CAUTION: 
these to modules should be always run in sequence, i.e., 

readpart [f n] 
decay [fnl] 

• prop [part] displays the properties of particle part 

• prrat [parti ,part2] gives the ratio of yields of particles parti and 
part 2 after the resonances have decayed 

• prratO [parti ,part2] gives the ratio of 'primordial' (before the de- 
cays) yields of particles parti and part2 

The following flags are used throughout the notebook to control the out- 
put and decays: 
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• listp (False or True) - when True, echos the properties of read parti- 
cles in the readpart module 

• describe (False or True) - controls the output of the readpart and 
decay modules. It is recommended to set it to True when running for 
the first time 

• debug (False or True) - when True, tests of the input files are performed 
while running the readpart and decay modules. In particular, a check 
is made if the decays proceed from the heaviest to the lightest particle 
(otherwise the DECAY OUT OF SEQUENCE warning is printed), the 
conservation of the baryon number and the electric charge is checked, 
and finally, the sum of the branching ratios is printed (of course, it 
should be equal to 1 within the assumed accuracy). It is recommended 
to set debug to True when running after modification of the input files. 

• yes3b (True or False) - when True, three-body decays are included. 
This is the normal mode. Switching off the three-body decays allows 
to test their importance 

• nomass (True or False) - when True, decays into product whose total 
mass is heavier than the mass of the decaying particle are allowed (!). 
There are a few decays like this in the Particle Data Tables, where the 
decaying particle is sufficiently wide. To switch off such decays, set 
nomass to False 

• threshmass (Infinity or number) - when Infinity, all decays are in- 
cluded, when set to a number, only particle lighter than threshmass 
are decayed 

• tracking (False or True) - when True, tracking of feeding to particle 
tracklab is made 

For convenience, particle labels negatives and positives are intro- 
duced. They include all negative or positive charge ground-state hadrons 
which are stable with respect to the strong interactions. 

Furthermore, we include the weak-decays feeding corrections as done by 
the STAR collaboration at RHIC. The STAR procedure excludes the feed- 
ing of pions from the decays of A(1115). We define labels of such pions as 
pi0139pluSTAR or pi0139minSTAR at the end of the decay module. Sim- 
ilar corrections can be defined for other ways of including the weak-decay 
corrections. 
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The format of the chi2.data file, containing the experimental informa- 
tion for the fitting, is 

parti part2 ratio stat.err syst.err 

More information on how to run Math SHARE is included in comments in 
the notebook itself. 

9 Conclusion 

We hope that SHARE will be useful both as an auxiliary tool in the analysis of 
experimental data and as a benchmark to verify more detailed ideas about the 
statistical approach to particle production in relativistic heavy-ion collisions. 
Also, the modules with the particle data and the resonance decays may be 
useful in their own right and interfaced to other models of particle production. 
We hope that given the versatility and modular structure of SHARE it will 
become a useful widely applied tool in study of hadron production in heavy- 
ion physics. 
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DATE : 20/07/04 TIME: 06.47.04 

UNITS : ENERGY IN GeV, DENSITY IN l/fm A 3 



temp 
lamq 
lams 
gamq 
gams 
lmi3 
lame 
game 
norm 
MU_B , MU_S 



.1373+/- 
1 .0591+/- 
1 . 02737723 

1 . 4294+/- 
1 .8735+/- 
1 .1014+/- 

1 . 
1 . 

2481 . 0208 



0.0102 limits: 
0.0439 limits: 



0.6479 limits: 
1.1862 limits: 
0.1457 limits: 



0.1000 
1 . 0000 

. 1000 
0.1000 
. 7500 



.1800 
10.0000 



2 . 1000 



0000 
2500 



0.0236606956 0.00417839215 



TOP 


BOTTOM 


THEORY 


EXP 


ERROR 


CHITERM 


Xil321min 


pi0139min 


. 


.0089 





.0093 


0. 


.0020 


-0 


.1954 


Xil321mnb 


pi0139min 


. 


.0074 





.0080 


0, 


. 0017 


-0. 


.3426 


Xil321mnb 


Xil321min 


. 


.8326 





.8530 


0. 


.1000 


-0 


.2043 


UM1672mnt 


Xil321mnt 


0. 


.1422 





.1438 


0. 


.0350 


-0 


. 0443 


Xil321min 


Lmlll5zer 


0, 


.1453 





.1870 


0. 


.4600 


-0. 


.0906 


Xil321mnb 


Lmlll5zrb 


0. 


.1572 





.2150 


0, 


. 5400 


-0 


. 1070 


Xil321min 


phl020zer 


0. 


.3162 





.3400 


0. 


.0700 


-0. 


.3397 


Lmlll5zrb 


Lmlll5zer 


0, 


.7696 





.7400 


0. 


.0500 





.5929 


Ka0492min 


pi0139min 


0. 


2073 





.1810 


0. 


.0250 


1 . 


. 0533 


Ka0492plu 


pi0139plu 


0. 


.1959 





.2010 


0. 


.0250 


-0 


.2053 


Ka0492plu 


pr0938plu 


2 . 


.0164 


2 


.4000 


0. 


.5000 


-0 


.7671 


Ka0492plu 


pr0938plb 


2 . 


. 9504 


3 


.0000 


0. 


.5000 


-0. 


.0992 


pr0938plb 


pr0938plu 


0. 


. 6834 





.7100 


0. 


.0600 


-0 


. 4426 


Ka0492plu 


Ka0492min 


1, 


.1185 


1 


.1030 


0. 


.0250 





. 6202 


Ka0492plu 


Ka0492min 


1, 


.1185 


1 


.1530 


0. 


.0400 


-0 


.8623 


Lml520zer 


Lmlll5zer 


0. 


.0421 





. 0220 


0. 


.0100 





.0000 


rho770zer 


pi0139min 


0. 


.0858 





.1000 


0. 


.0200 


-0 


.7076 


f 00980zer 


pi0139min 


. 


.0065 





.0200 


0. 


.0100 


-1 . 


.3547 


Ka0892zer 


Ka0492min 


. 


.2424 





.2400 


0. 


.0100 





. 2413 


phl020zer 


Ka0892min 


0. 


.5803 





.5950 


0. 


.2400 


-0 


.0611 


D11232plp 


pr0938plu 


. 


.2002 





.2400 


0. 


.0370 


-1 . 


. 0762 



DATA POINTS: 21 PARAMETERS: 5 DoF : 16 
TOTAL CHI**2/DEG. OF FREEDOM 0.457254927 
SIGNIFICANCE OF FIT: 0.966676891 
[ (S-SBAR) / (S+SBAR) ] -5 . 5 6 8 7 8 37E-1 7 



Figure 4: A typical fitratios.out file 
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